NASA Technical Memorandum 83548 
AIAA-84-0531 


Computation of Viscous Flow in 
Curved Ducts and Comparison With 
cLk* erimeiital Data 


t MAS A Vja- 83 54 8) COSPOTATJCN Of VISCOJS FlOU N84- 13494 

IN CO V£D DUCIS A 80 CCflPAfilSON WITH 

E1BEB H2NIAL DATA (NASA) 22 p HC A02/HP AO 1 

CSCh 20D Unclas 

G3/34 42739 



Charles E. Towne 
Lewis Research Center 
Cleveland , Ohio 


Prepared for the 

Twenty-second Aerospace Sciences Meeting 
sponsored by the American Institute of 
Aeronautics and Astronautics 
Reno, Nevada, January 9-12, 1984 





E— 1922 


COMPUTATION OF VISCOUS FLOW IN CURVED DUCTS AND COMPARISON 
WITH EXPERIMENTAL DATA 
Charles E. Towne 

National Aeronautics and Space Administration 
Lewis Researci Center 
Cleveland, Onio 44135 


SUMMARY 

A three-dimensional analysis for fully viscous subsonic internal flow is 
eva’uated. The analysis, designated PEPSIG, solves an approximate form of the 
Navier-Stokes equations by an implicit spatial marching procedure. Results of 
calculations are presented for laminar flow through two different circular 
cross-sectioned 180 degree bends, and for laminar and turbulent flow through 
circular and square cross-sectioned 22.5-22.5 degree S-ducts. Quantitative 
comparisons with experimental data are shown for all cases. Special emphasis 
is placed on verifying the ability of the analysis to accurately predict the 
distorted flow fields resulting from pressure-driven secondary flows. 


INTRODUCTION 

The presence of three-dimensional subsonic internal flow is relatively 
common in many of today's airbreathinq propulsion systems. Curved centerlines 
and changes in cross-sectional shape are often present, resulting in the gen- 
eration of complex three-dimensional secondary flows and significant flow dis- 
tortion. In addition, viscous effects are usually important, with the boundary 
layer thickness often comparable to the duct radius. These phenomena have im- 
portant effects on overall propulsici system performance. 

Conventional boundary layer methods will not work for these complex flows. 
A full Navier-Stokes solution could be used, but would require large amounts 
of computer time and storage. However, by making certain approximations to the 
Navier-Stokes equations, such as neglecting streamwise diffusion, a set of 
equations can be derived for fully viscous internal flow that can be solved by 
forward marching in space. An analysis method using tnese equations to compute 
three-dimensional subsonic flow through curved ducts with super-elliptic cross- 
sections was developed by Briley and McDonald (ref. 1) and by Levy, McDonald, 
Briley, and Kreskovsky (ref. 2), and was designated PEPSIG. A preliminary 
evaluation of the PEPSIG method was presented by Towne and Anderson (ref. 3). 
However, in that evaluation only qualitative comparisons were made with data. 

In addition, ieveral changes have since been made to the analysis by Levy, 
Briley, and Mc r )o"ald (ref. 4). 

The objective of the present study, therefore, was to evaluate and verify 
the most recent version of the PEPSIG analysis by making detailed quantitative 
comparisons between predicted results and benchmark experimental data. Special 
emphasis was placed on testing the ability of the analysis to accurately pre- 
dict the distorted flow field resulting from pressure-driven secondary flows. 



GOVERNING EQUATIONS 


In this analysis, the flow is computed by a single sweep spatial marching 
procedure which solves an approximate form of the Navier-Stokes equations. It 
is assumed that the flow is primarily in the direction of the duct centerline, 
with transverse secondary flow. This allows two basic assumptions to be made. 
The first is that second derivatives in the primary flow direction are negli- 
gible. The second is that the pressure in the primary, or streamwise, momentum 
equation can be represented by the sum of a known three-dimensional pressure 
field and a one-dimensional correction computed as part of the marching proce- 
dure to account for viscous blockage. When these two assumptions are applied 
to the Navier-Stokes equations, a set of equations can be derived that can be 
solved by forward marching in the primary flow direction. 

The details of the derivation of the governing equations have been presen- 
ted elsewhere (refs. 1, 2, and 4). For completeness, however, a brief discus- 
sion is presented here. The equations are first written in a centerline-based 
orthogonal reference coordinate system denoted by xj, x^, and X 3 , with cor- 
responding velocities v s , w s , and Up. The X 3 direction is the primary flow 
direction, and is aligned with the duct centerline at each marching step. The 
x} and x^ directions are normal to the duct centerline and define the trans- 
verse, or secondary, flow plane. The primary momentum equation is then given 
by 
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This equation, and the others to be presented in this section,, are taken 
directly from reference 4. Here h is the metric scale coefficient associated 
with the primary direction X 3 . (The metric scale coefficients for the trans- 
verse directions are both unity). Also, Pj = Pi (xj.xp.xj) is the known three- 
dimensional pressure field and p v * P v (x 3 ) is the unknown one-dimensional 
viscous blockage correction. In this study, p,- was computed using a three- 
dimensional potential flow analysis, thus bringing elliptic effects due to 
geometry into the solution. 

The secondary velocities v $ and w $ are split into irrotational 
and rotational components. The continuity equation can then be written as 
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Here q> is a scalar potential corresponding to the irrotational component of 
the secondary velocity. 


The transverse momentum equations are cross-differentiated and combined to 
form an equation for the transport of streamwise vorticity in the primary flow 
direction, given by 
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Here 9 is the stream function, or vector potential, corresponding to the 
rotational component of the secondary velocity, and a is the vorticity in 
the X3 direction. 


The fourth equation is simply the definition of a written in terms of 9. 



+ a = 0 


The energy equation is eliminated by assuming constant total enthalpy. For tur- 
bulent flow, an algebraic mixing length turbulence model is used. 

The above equations are transformed into a non-orthogonal body-fitted co- 
ordinate system. They are then solved by forward marching in the streamwise 
direction from an initial station, where the flow profiles are known, by an im- 
plicit finite-difference technique. The boundary conditions that are applied 
result in no-slip at the walls, and symmetry at the symmetry plane. Details of 
the solution procedure can be found in references 1, 2 , and 4. 


COMPUTED RESULTS 

The results of the PEPSIG analysis have been compared with experimental 
data for a variety of geometric configurations and flow conditions. In this 
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section, results are presented and discussed for the following cases: (1) in- 

compressible laminar flow through two different circular cross-sectioned 180 
degree bends tested by Agrawal, Talbot., and Gong (ref. 5); (2) incompressible 
laminar and turbulent flow through a circular cross-sectioned 22.5-22.5 S-duct 
tested by Taylor, Whitelaw, and Yianneskis (ref. 6); and (3) incompressible 
laminar and turbulent flow through a square cross-sectioned 22.5-22.5 S-duct 
also tested by Taylor, Whitelaw, and Yiannes is (ref. 7). All of these data 
are of benchmark quality, with streamwise and cross-flow velocities measured 
using laser-doppler velocimetry. PEPSIG results for other configurations are 
presented by Anderson (ref. 8) and by Vakili, et al . (ref. 9). 


180 Degree Bend 

The two geometries tested in reference 5 are shown in figure 1. Both have 
a cross-section radius r of 1.905 cm (0.75 in.). They had centerline radii 
R of 13.335 cm (5.25 in.) and 38.1 cm (15.0 in.), corresponding to radius 
ratios R/r of 7 and 20, respectively. 

R/r = 7 . - The first results to be presented are for the R/r = 7 bend. 
The flow was laminar, with a Reynolds number _Re of 242 based on the cross- 
section radius r and the average velocity u, corresp ondi ng to a Dean number 
De of 183. The Dean number is defined as De = 2 Re y/r/R. With a 20x20 mesh 
in half the cross-section (taking advantage of symmetry), and 181 streamwise 
marching steps, 1.7 minutes of CPU time were needed on a Cray-1 computer. 

The computed secondary flow pattern after 60 degrees of turning is shown 
in figure 2. Each arrow in this figure represents the direction and magnitude 
of the secondary velocity at the tail of the arrow. This is the classic pat- 
tern for flow in a curved pipe. The low energy flow near the wall migrates 
circumferentially away from the outside, or pressure side, of the bend and 
toward the inside, or suction side. The higher speed flow in the core region 
moves toward the outside of the bend due to centrifugal effects. A pair of 
counter-rotating vortices is thus established. This secondary flow pattern 
persists through the length of the duct. 

The resulting flow distortion is shown in figure 3 in the form of contours 
of constant streamwise velocity. These results are presented at the same sta- 
tions used in the experiment. The boundary layer initially thickens slightly 
on the outside of the bend due to the pressure increase there as the flow en- 
ters the bend. By the second station, at 30 degrees, the effect of the second- 
ary flow is apparent in the thickening of the boundary layer on the inside of 
the bend. By 60 degrees of turning the secondary flow has caused significant 
flow distortion. A pocket of low speed flow has formed on the inside of the 
bend, displacing the region of maximum velocity toward the outside of the bend. 
The results at the final two stations are similar. 

In the experiment of reference 5, data were taken at five streamwise sta- 
tions. At each station, streamwise velocity profiles were measured along the 
five data lines shown in the inset to figure 4. The computed streamwise ve- 
locity profiles along these data lines are compared with the data in figure 4. 
Computed results are shown for both 20x20 and 40x40 cross-sectional meshes; 
however, the effect of greater mesh resolution for this case is negligible. 
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For both meshes, the agreement between the PEPSIG results and the data is 
excellent. 

The computed and experimental secondary velocity profiles are presented 
in figure 5. These results, both computed and experimental, are for slightly 
different flow conditions. The Reynolds and Dean numbers were 183 and 138, 
respectively. The distance between adjacent data lines in figure 5 corresponds 
to a v/u value of 0.36. The agreement between theory and experiment is gen- 
erally very good. 

R/r = 20 . - For the R/r = 20 case, the flow was laminar with a Reynolds 
number of 1263 and Dean number of 565. Several different meshes were used. 

For a 50x50x226 mesh, 13.1 minutes of CPU time were used on a Cray-1 computer. 

The computed secondary flow pattern after 84 degrees of turning is shown 
in figure 6. The basic pattern is similar to that seen in the R/r = 7 case, 
but the boundary layers are thinner, the secondary velocities are not as large, 
and the vortices are tighter. The computed streamwise velocity contours are 
shown in figure 7. Again, these results are qualitatively similar to those in 
the R/r = 7 case. However, the flow is more distorted, and becomes distorted 
after less turning. 

In figures 8(a) and (b), the computed streamwise velocity profiles are 
compared with the data. (The results along data line 5 are not shown because 
of problems encountered in interpolating the computed results in the high ve- 
locity gradient region near the wall.) In figure 8(a), results for both 25x25 
and 50x50 cross-sectional meshes are shown. A total of 226 points were used in 
the streamwise direction, with ae = 0.5 degrees for e between 0 and 45 deg- 
rees, and ae =1.0 degree for e greater than 45 degrees. The agreement is 
generally very good for both meshes. In fact, the 25x25 mesh gives better re- 
sults than the 50x50 mesh at e = 108 degrees. However, the coarser me' . is 
unable to resolve the velocity peak near the inner wall at e = 35 degrees seen 
along data lines 2 and 3. In figure 8(b), the results for the same 50x50x226 
mesh are shown along with results for a 50x50x181 mesh with ae = 1.0 degree 
for the entire duct. With the larger streamwise step size, the analysis is 
again unable to resolve the inner wall velocity peak at e = 35 degrees. 

These results indicate that proper mesh resolution in both the transverse and 
streamwise directions is critical in correctly predicting the initial develop- 
ment of the secondary flows. 

Finally, in figure 9 the computed and experimental streamwise velocity 
contours are compared after 84 degrees of turning. The computed results are 
for the 50x50x226 mesh. The agreement between theory and experiment is 
excellent. 


Circular S-Duct 

The circular cross-sectioned S-duct geometry tested in reference 6 con- 
sisted of two 22.5 degree circular arc bends, as shown in figure 10. The cross- 
section diameter D was 48 mm and the centerline radius R was 336 mm, for a 
ratio R/D of 7. 
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Laminar flow . - For the laminar flow case, _the Reynolds number based on 
cross-section diameter D and average velocity u was 790, corresponding to 
a Dean number of 211. With a 20x20x80 mesh, this case required 0.8 minutes of 
CPU time on a Cray-1 computer. 

The computed secondary velocities near the inflection and exit planes are 
shown in figure 11. In the first bend the pressure-driven secondary flow pat- 
tern typical of curved pipes is set up. In the second bend, the cross-flow 
pressure gradients change sign, tending to attenuate the secondary flows set up 
in the first bend. In fact, by the time the exit plane is reached, a vortex- 
type motion has been set up in the lower half of the duct (i.e., near the inner 
wall of the second bend) that is in the opposite sense as that in the first 
bend. However, as noted in references 6, 7, and 10, the distortion of the 
streamwise flow in the first bend causes the generation of streamwise vorticity 
in the top half of the second bend that is in the same sense as that already 
present, in accordance with the Squire and Winter formula for the development 
of streamwise vorticity. The vortices set up in the first bend are thus sus- 
tained, or even intensified, in the second bend and become concentrated near 
the top of the duct. 

The effect of the secondary flows on the streamwise flow is shown by the 
streamwise velocity contours in figure 12. The boundary layer along the bottom 
of the duct initially thickens due to the adverse pressure gradient there as 
the flow enters the duct. By the time the second station is reached, just 
ahead of the inflection plane, the secondary flows have caused the boundary 
layer in the top half of the duct to thicken significantly. Through the in- 
flection region an adverse pressure gradient exists along the top of the duct, 
and the analysis predicts a small region of separated flow. This separation 
region has no appreciable effect on the rest of the flow, and the "FLARE" ap- 
proximation (ref. 11) is used to allow the streamwise marching analysis to con- 
tinue. In the second half of the duct the concentration of the vortices near 
the top causes a severely distorted flow field. In addition, the secondary ve- 
locities near the wall in the bottom half of the duct cause the boundary layer 
there to thicken. The flow at the duct exit is as distorted, or more distor- 
ted, than it would be in a unidirectional bend with the same total amount of 
turning. 

In figure 13, the computed streamwise velocity profiles in the symmetry 
plane are compared with the data. Analytical results are shown for both 
20x20x80 and 40x40x80 meshes. Both give very good agreement with the data, 
and the effect of greater mesh resolution is negligible. 

Turbulent flow . - The Reynolds number for turbulent flow in the circular 
S-duct was 48,000, corresponding to a Dean number of 12,828. With a 50x50x80 
mesh, this case required 4.8 minutes of CPU time on a Cray-1 computer. 

The computed secondary velocities near the inflection and exit planes are 
shown in figure 14. The physical phenomena here are the same as in the laminar 
flow case. However, the boundary layers are much thinner, resulting in gener- 
ally lower secondary velocities. The vortices are also somewhat more concen- 
trated and nearer the top of the duct. The streamwise velocity contours are 
shown in figure 15. Again, except for the thinner boundary layers, these re- 
sults are qualitatively similar to the laminar results. 
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In figure 16, the computed and experimental streamwise velocity profiles 
in the symmetry plane are compared with «.he data. In contrast to the laminar 
case, mesh resolution for the turbulent case has a significant effect on the 
results. With the 50x50x80 mesh, the PEPSIG results agree very well with the 
data. With the 25x25x80 mesh, however, the thin boundary layers and secondary 
flow vortices are not resolved well enough for the analysis to predict the 
distorted flow region at the top of the duct. 


Square S-Duct 

The square cross-sectioned S-duct geometry of reference 7 also consisted 
of two 22.5 degree circular arc bends, as shown in figure 10. The cross-section 
width 0 was 40 mm and the centerline radius R was 280 mm, for a ratio R/D 
of 7, the same as the circular cross-sectioned S-bend. Since the PEPSIG com- 
puter code is presently limited to super-elliptic cross-sections, the square 
cross-section was represented by a super-ellipse with an exponent of 10. 

Laminar flow . - For the laminar flow case, the Reynolds number based on 
cross-section width 0 and average velocity u was 790, corresponding to a 
Dean number of 211, the same as in the circular cross-sectioned laminar flow 
case. With a 20x20x111 mesh, this case required 1.2 minutes of CPU time on a 
Cray-1 computer. 

The computed secondary velocities near the inflection and exit planes are 
shown in figure 17. The important physical phenomena for this case are the 
same as for the circular cross-sectioned case, and the secondary flow patterns 
are very similar. A pair of counter-rotating vortices is set up in the first 
bend. In the second bend, these vortices become concentrated near the top of 
tne duct, and the near-wall secondary velocities in the lower half of the duct 
change direction. 

The computed streamwise velocity contours are shown in figure 18. These 
results are, again, very similar to those in the circular cross-sectioned 
S-duct. In the first bend, the boundary layer in the top half of the duct 
thickens because of the secondary flows. In the second bend, the concentration 
of the vortices near the top of the duct again causes a severely distorted flow 
field. Small regions of separated flow were predicted for this case along the 
top of the duct near the inflection plane and along the bottom of the duct near 
the exit plane, where adverse pressure gradients exist. 

The computed and experimental streamwise velocity profiles in the symmetry 
plane are compared in figure 19. Both the 20x20x111 and 40x40x111 mesh results 
agree very well with the data. The 20x20x111 case shows some undesirable wig- 
gles, however, at the last two stations in the middle of the duct where the 
mesh resolution is poor. 

Turbulent flow . - The Reynolds and Dean numbers for the turbulent flow 
case were 40,000 and 10,690, respectively. With a 50x50x111 mesh, this case 
required 7.7 minutes of CPU time on a Cra.y-1 computer. 

The computed secondary velocities ear the inflection and exit planes are 
shown in figure 20. As in the circular cross-sectioned case, these results are 
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basically similar to the corresponding laminar results, but with thinner bound- 
ary layers and generally lower secondary velocities. However, the vortices at 
the top of the duct have essentially disappeared by the time the exit plane is 
reached. In addition, in the second bend the streamwise velocity contours for 
turbulent flow, shown in figure 21, are shaped differently than for laminar 
flow. Near the top of the duct they have a gentle "W" shape, with the low 
speed flow extending furthest into the core region a short distance away from 
the symmetry plane. In the laminar flow case, in contrast, the low speed flow 
extended furthest into the core region exactly in the symmetry plane. The 
shape of these contour lines is in agreement with that found experimentally. 

In figure 22 the computed and experimental streamwise velocity profiles 
in the symmetry plane are compared with the data. Computed results are shown 
for both 25x25x111 and 50x50x111 meshes. The 50x50x111 mesh gives slightly 
better agreement with the data than the 25x25x111 mesh, although the differen- 
ces are small. 


CONCLUDING REMARKS 

For the configurations studied, the PEPSIG analysis gave quantitatively 
accurate results for both laminar and turbulent flow when sufficient mesh reso- 
lution was used. In particular, the analysis correctly predicted the complex 
physical phenomena and distorted flow fields that occur in S-shaped ducts. For 
the cases with thick inlet boundary layers, 20x20 cross-sectional meshes were 
sufficient. However, for the thin boundary layer cases a 50x50 cross-sectional 
mesh was generally needed for accurate results. Proper mesh resolution in both 
the transverse and streamwise directions was especially important in regions 
where the secondary flows were initially developing. 

The PEPSIG analysis is very fast compared to the alternative, a full 
Navier-Stokes calculation. For the results presented in this paper, the cur- 
rent PEPSIG computer code processed 600-700 grid points per second on the 
Cray-1 computer. (The analysis is approximately 10 times slower on an IBM 370/ 
3033). In reference 12, a time o f 0.028 seconds per grid pent per time step 
is reported for a three-dimensional problem solved by an efficient full Navier- 
Stokes procedure on an IBM 370/3033. Assuming 100 time steps for a converged 
solution, and a speed-up factor of 10 between the IBM and Cray, the Navier- 
Stokes procedure would process 3.6 grid points per second. The PEPSIG analysis 
is thus approximately 150-200 times faster, for comparable accuracy, than an 
efficient Navier-Stokes procedure. 
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Figure 4 - Computed and experimental streamwise velocity profiles in R/r * 7 duct 
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Figure 5. - Computed and experimental cross-flow velocity profiles in 
R / r ■ 7 duct 
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Figure 6. - Computed secondary flow at 0 « 8^ 
in R/r * 20 duct. 



Figure 7. - Computed streamwise velocity contours in K/r - 20 duct. 
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Figure 9. - Computed and experimental streamwise velocity contours at 8 - 84P in R/r - 20 
duct 



Figure 10. - S -duct configuration. 
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Figure U. - Computed secondary velocities for laminar flow in 22. 5-22. 5 circular S-duct 



Figure 12. - Computed streamwise velxity contours 
for laminar flow in 22. 5 - 22. 5 circular S-duct 
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Figure 13. - Computed and experimental streamwise velocity profiles in symmetry 
plane for laminar flow in 22. 5 - 22. 5 circular S-duct 
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Figure 14 * Computed secondary velocities for turbulent flow in 22. 5 - 22, 5 circular S-duct 
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Figure 15. - Computed streamwise veixity contours 
for turbulent flow in 22. 5 - 22. Scircular S-duct 
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Figure 16. - Computed and experimental 
streamwise veixity profiles in symmetry 
plane for turbulent flow in 22. 5 - 22. 5 
circular S-duct 
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Figure 17. - Computed secondary velocities for laminar flow in 22 . 5 - 22. 5 square S-duct 



Figure 18. - Computed streamwise velocity contours for laminar flow in 
22.5 - 22.5 square S-duct. 
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Figure 19. - Computed and experimental streamwise ve- 
locity profiles in symmetry plane for laminar flow in 
22. 5 - 22. 5 square S-duct 
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Figure 20. - Computed secondary velocities for turbulent flow in 22. 5 - 22, 5 square S-duct 
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Figure 21, - Computed streamwise velocity contours for turbulent flow in 
22, 5 - 22. 5 square S -duct 
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Figure 22. - Computed and experimental stream- 
wise velxity profiles in symmetry plane for 
turbulent flow in 22. 5 - 22. 5 square S-duct 




